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NOTICE 


When  U.S.  Government  drawings,  specifications,  or  other  data  are  used  for 
any  purpose  other  than  a  definitely  related  government  procurement  operation, 
the  government  thereby  incurs  no  responsibility  nor  any  obligation  whatsoever, 
and  the  fact  that  the  government  may  have  formulated,  furnished,  or  in  any  way 
supplied  the  said  drawings,  specifications,  or  other  data,  is  not  to  be 
regarded  by  implication  or  otherwise,  or  conveying  any  rights  or  permission  to 
manufacture,  use,  or  sell  any  patented  invention  that  may  in  any  way  be 
related  thereto. 
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This  final  report,  A  Diffuser  Heat  Transfer  Code,  was  prepared  by  Dr  Gale 
H.  Buzzard  of  Duke  University  during  the  summer  of  fiscal  year  1983  as  an  in- 
house  project  while  serving  at  the  Air  Force  Rocket  Propulsion  Laboratory, 
Edwards  Air  Force  Base,  CA,  as  a  visiting  professor.  The  project  manager  for 
the  Air  Force  Rocket  Propulsion  Laboratory  was  Dr  Philip  Kessel. 
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1 .  INTRODUCTION 


A  computer  program  has  been  developed  that  performs  a  thermal  analysis  of 
a  water  jacketed  rocket  motor  test  diffuser.  The  program  has  been  developed 
to  handle  the  requirements  of  the  particle  laden  plusr  associated  with  a 

metallized  solid  propellant  but  is  also  capable  of  handl’og  a  particle  free 

l  2 

plume.  The  program  combines  the  earlier  work  of  Trout  and  Me  Cay  ,  Pergament  , 
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and  Kessel'.  The  end  result  is  a  Diffuser  Heat  Transfer  code  (D1ST)  which 
corrects  several  of  the  shortcomings  of  the  Rocket  Engine  Thermal  Analysis 
Program  (REDTAP)  created  by  Trout  and  McCay  and  includes  severa'  areas  not 
treated  by  the  earlier  code.  Included  among  these  areas  are  .adial 
temperature  gradient  within  the  diffusei  wall,  impingement  heat  transfer  at 
the  point  of  plume  attachment,  an  improved  model  for  the  particle  impingement 
accommodation  coefficient,  and  particle  debris  shielding. 

In  addition  to  the  development  cf  the  !)u>”  code,  the  77-inch  diffuser 
located  at  the  Air  Force  Rocket  Propulsion  Laboratory  (AFRPL1  Test  Area  1-42 
was  instrumented  to  record  water  side  wall  temperature  and  coolant  temperature 
at  selected  sites  along  the  initial  seven  feet  of  the  diffuser  during  routine 
test  firings.  It  was  anticipated  that  data  would  be  available  early  in  the 
project.  The  principal  thrust  of  the  present  work  was  to  have  been  the 
writing  of  the  DHT  code  followed  by  a  comparison  of  the  DHT  predictions  with 
the  experimental  data  gathered  at  Test  Area  1-42.  It  was  aniicipated  that  the 
code  would  evolve  as  the  body  of  experimental  data  grew  and  that  the  end 
result  would  be  a  diffuser  heat  transfer  code  supported  by  a  bo.lv  of 
experimental  data.  Unfortunately,  the  test  firings  to  date  have  involved 
motors  that  iiave  been  too  small  and/or  burn  times  that  have  been  too  short  to 
provide  data  tiiat  ere  useful  in  validating  the  DHT  code.  The  quality  of  the 
data  thus  far  has  been  excellent  and  there  Is  the  promise  of  useful  data  at  a 
future  date,  but  much  remains  to  be  done  in  terms  of  validating  the 
predictions  of  DHT. 

As  was  the  case  with  the  earlier  code,  DHT  relies  on  the  AFRPL  Solid 

4 

Performance  Program  (SPP)  and  a  modification  of  the  Joint  Army  Navy  NASA  Air 
Force  UANNAF)  Standardized  Plume  Flow  Field  Model  (SCIPPY;  ’  to  provide  the 
flow  field  data  within  the  diffuser;  DHT  incorporates  the  Inter-Agency 
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Chemical  Rocket  Propulsion  Croup  Turbulent  Boundary  Layer  code  (T8L)^  as  a 
subroutine  to  handle  the  gas  side  convection  heat  transfer.  Unfortunately, 

the  dependency  of  DHT  upon  SCIPPY  places  a  current  restriction  on  the 
usefulness  of  DHT.  SCIPPY  is  still  in  che  developmental  stage  and  currently 

is  unable  to  handle  rocket  motor/diffuser  size  combinations  that  involve  large 
disparities  between  the  exit  diameter  of  the  motor  and  the  inlet  diameter  of 
the  diffuser.  Hopefully,  as  SCIPPY  evolves  toward  its  final  form,  these 
problems  wil)  disappear. 
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1.  STATEMENT  OF  THE  PROBLEM 


2.1  GENERAL 


The  computer  code  developed  and  all  discussions  that  follow  deal 
specifically  with  the  77-inch  diffuser  located  in  the  AFRPL  Test  Area  1-42. 
The  code,  however,  is  written  to  perform  a  thermal  analysis  for  any  similar 
diffuser  and  could  he  easily  modified  to  handle  most  water  jacketed 
conf igurat ions , 
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Shown  in  Figure  1  is  a  schematic  of  the  diffuser.  The  diffuser  has  a 
uniform  diameter  inlet  section  followed  by  a  conical  transition  to  a  second 
uniform  diameter  section  and  conical  expansion.  This  final  expansion  connects 
to  a  plenum  which  removes  the  rocket  motor  exhaust  gases  and  maintains  the 
reduced  pressure  necessary  to  simulate  altitude  conditions. 

The  diffuser  is  fabricated  from  ASME-SA-285-C  steel  and  has  a  0.50-inch 
inner  wall  which  forms  the  containment  for  the  exhaust  plume.  The  water 
jacket  is  formed  by  this  and  a  0.375-inch  outer  wall.  These  two  walls  are 
separated  by  2.5  x  ?  j  x  0.5-inch  angle  wound  with  a  5.75-inch  pitch  quadruple 
lead  that  results  in  four  parallel  coolant  passages  approximately  5.25  x  2.75 
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inches.  These  angle  members  are  welded  to  the  inner  wall.  No  attempt  has 
been  made  to  analyze  the  thermal  path  added  by  these  angle  members.  One  can 
assume  that  they  will  provide  additional  cooling  of  the  inner  wall  out  the 
extent  of  this  effect  is  indeterminate.  There  is  the  added  complication  of  a 
nominal  0.25~inch  radial  clearance  between  the  inner  assembly  and  the  outer 
wall.  Since  the  inner  assembly  floats  within  the  outer  wall,  the  resulting 
radial  clearance  can  range  from  0  to  0.50  inches.  It  is  assumed  that  this 

condition  does  not  short  circuit  the  helical  path  of  the  water  jacket.  This 
is  a  question  that  must,  be  addressed  as  experimental  data  become  available. 

The  heat  load  on  the  diffuser  is  comprised  of  the  convective  load  from 

the  exhaust  gases  plus  the  various  particle  related  heat  fluxes.  The  particles 
ca^ry  with  them  a  very  significant  quantity  of  thermal  energy  as  a  result  of 
their  heat  capacity  and  elevated  temperature  in  reference  to  the  gas  side  wall 
tempei store  of  the  diffuset.  They  also  carry  a  very  significant  quantity  of 
kinetic  energy.  If  as  they  impact  the  wall,  an  appreciable  portion  of  either 

of  these  energies  is  transferred  to  the  diffuser  wall,  a  very  severe  heat  load 

will  result.  Crucial  to  a  valid  diffuser  model  is  the  selection  or 
development  of  a  particle  impingement  model  that  adequately  handles  the 
exchange  of  these  two  forms  of  energy.  Radiant  exchange  from  the  particles  to 
the  wall  is  a  non-negligible  but  distinctly  second  order  heat  load. 

The  diffuber  wall  must  obey  the  unsteady  heat  equation.  It  is  convenient 
to  note  that  the  wall  is  thin  compared  with  the  diffuser  radius  and  to  write 
the  he.it  equa* ion  in  two-dimensional  Cartesian  coordinates.  As  such,  the  wall 
temperature  is  governed  by 
2.  2,, 

J_J  +  J..J  =  (0c/k)  6T/6t  (1) 

c  2  i  2 
x  y 

where  x  and  y  are  measured  parallel  and  normal  to  the  diffuser  wall  and  i>  ,  C 
and  k  are  the  density,  specific  heat  and  thermal  conductivity  of  the  diffuser 
wall. 

The  fluid  Clow  within  the  water  jacket  is  assumed  to  be  one-dimensional 
constant  property  steady  flow.  It  is  assumed  th*,t  the  flow  rate  is  known  and, 
therefore,  the  local  fluid  velocity  is  a  simple  function  of  the  local  water 
jacket  cross-sectional  area.  The  temperature  distribution  within  Che  water 
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jacket  is  assumed  to  be  a  one-diraens ional  axial  transient  superimposed  upon 
the  steady  flow  and  involving  negligible  axial  conduction.  As  such,  the  water 
jacket  temperature  is  governed  by 


u  ~  +  (2  !T  Rh/  p  CA)  (T  -  T  .  .  )  =  -  6T/6t  (2) 

^  X  Wfl  L  1 

where  U  is  the  axial  velocity  of  the  coolant,  R  is  the  outer  radius  of  the 
inner  wall,  r  and  C  are  the  density  and  specific  heat  of  the  coolant,  A  is 
the  axial  cross-sectional  area  of  the  water  jacket,  and  h  is  the  water  side 
fill),  coefficient.  The  outer  wall  of  the  water  jacket  is  treated  as  an 
adiabatic  surface.  Equations  1  and  2  are  solved  using  explicit  finite 
difference  techniques. 


2.2.  PARTICLE  IMPINGEMENT 

The  aluminum  oxide  particles  contained  within  the  exhaust  of  a  metallized 
solid  propellant  rocket  motor  carry  with  them  a  considerable  quantity  of 
thermal  and  kinetic  energy.  It  is  convenient  to  measure  the  thermal  energy 
relative  to  the  gas  side  wall  temperature  of  the  diffuser  and  to  partition  the 
kinetic  energy  into  a  component  resulting  from  the  velocity  parallel  to  the 
diffuser  wall  and  a  component  resulting  from  the  velocity  normal  to  the  wall. 
In  this  form  the  potential  heat  load  caused  by  the  particles  impinging  on  the 
diffuser  wall  may  be  represented  by 


m  C  (T  -  T  )  +  in  U2/2  +  m  V2/2 
ppp  w  pp  pp 


where  m  is  the  mass  Cl  iw  of  particles  impinging  upon  the  wall,  is  tlu* 

specific  heat  of  the  particles,  T  is  the  temperature  of  the  particle,  T^  is 
the  gas  side  wall  temperature  of  the  diffuser,  U  is  the  velocity  of  the 
particle  parallel  to  the  wall  and  is  the  velocity  of  the  perticle  normal  to 
the  wall. 

I r  is  common  practice  to  quantify  the  particle/wal i  interaction  in  terms 
of  three  accommodation  coefficients  (C^,,  and  Cy)  which  define  the  fraction 
of  each  energy  component  that  is  transferred  to  the  diffuser  wall. 
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Introducing  this  concept,  the  heat  load  on  the  diffuser  oall  due  to  particle 
impingement  is  „iven  by  the  following  expression 


=  m  C_C  (1 
P  [  T  P 


T  -  T  )  +  C..U^/2  +  C„V2/2 

P  w  U  p  V  p 


Evaluation  of  these  accommodation  coefficients  is  in  large  measure  a 
question  of  particle  behavior  upon  impact  with  the  wall.  If  the  particles 
adhere  to  the  wall,  all  three  accommodation  coefficients  are  unity  and  the 
impingement  heat  flux  will  be  the  dominant  heat  load  on  the  diffuser.  Such  an 
assumption  would  be  a  very  safe  estimate  of  the  maximum  heat  flux  but  if 
overly  cautious  would  preclude  the  testing  of  rocket  motors  that  could  in 
reality  be  safely  tested  within  the  facility.  Particle  impingement  with  the 
wall  can  be  expected  to  occur  at  a  relatively  shallow  angle.  This  lends 
credence  to  an  assumption  that  the  particles  do  not  adhere  to  the  wall  and 
that  the  thermal  accommodation  coefficient  ;s  close  to  zero.  A  further 
consequence  of  this  assumption  would  be  that  the  momentum  of  the  particle 
parallel  to  the  surface,  and  therefore  that  component  of  the  kinetic  energy, 
will  be  conserved.  This  would  lead  to  a  C.j  equal  to  zero.  Visual  inspection 
following  two  Super  BATES  firings  revealed  no  evidence  of  significant  particle 
deposition  on  the  diffuser  wall. 

The  transfer  of  the  component  of  kinetic  energy  normal  to  the  surface  from 

the  particle  to  the  surface  can  be  related  to  the  coefficient  of  restitution 

.  .  .  3 

for  the  collision.  After  a  review  of  the  limited  data  available,  Kesael 

suggests  the  use  of  a  coefficient  of  restitution  equal  to  (1  -  B/90),  where  B 

is  the  angle  of  impact,  as  measured  in  degrees,  between  the  velocity  vector 

and  the  normal  to  the  surface.  This  leads  to  the  following  expressions  where 


V  is  the  component  of  the  velocity 
denotes  conditions  following  impact. 


normal 


to  the  surface  and  the  prime 


V  /V  =  1  -  B/90 
P  P 


<KE'/KE)non»al  *  (1  *  B/90>2 
(  AKE/KE)normal  =  (B/90){2  -  B/90) 


The  decrease  in  the  normal  compone;  t  of  the  kinetic  energy  places  an  upper 

hound  on  tiie  energy  transferred  to  the  surface.  Unless  the  particle  adheres 
to  the  hurl.tce,  a  portion  of  this  energy  will  be  carried  away  as  an  increase 

m  the  internal  energy  of  the  particle.  Citing  limited  data  that  support  an 

accommodation  coefficient  of  0.35  to  0.70  for  normal  impact  and  noting  that 
the  quantity 

(B/90)(2  -  B/ 90) 

is  approximated  within  t  IX  by  1.15  SIN  B  for  B  less  than  or  equal  to  40 

degrees,  Kessel  suggests  the  use  of  an  accommodation  coefficient 

Cy  =  0.8  SIN  B  (4) 

The  heat  load  associated  with  particle  impingement  is  handled  within  DHT 
as  per  Equation  3  with  the  user  allowed  to  specify  any  desired  set  of  values 
for  the  accommodation  coefficients.  SCIPPY  will  provide  DHT  with  local  values 
for  B  and  an  option  is  provided  that  allows  the  use  of  Equation  4  along  with 
the  ability  to  scale  the  coefficient  of  SIN  B  up  or  down  at  will.  In 
reporting  on  data  gathered  at  Arnold  Engineering  Development  Center  (AEDC)  on 

3 

an  instrumented  diffuser,  Kessel  cites  modest  agreement  between  the 
experimental  data  and  the  predictions  of  REDTAP  using  accommodation 
coefficients  of  0,  0,  and  0.8  SIN  B  along  wit!  a  specified  average  value  of  B 
equal  to  22  degrees. 

2.3  PARTICLE  RADIATION 

No  attempt  's  made  in  the  present  work  to  alter  the  radiation  model 
developed  within  REDTAP  by  Trout  and  McCay.  This  is  a  very  simplistic  model 
that  places  a  believable  upper  limit  on  the  contribution  of  particle  radiation 
to  the  heat  load  on  the  diffuser.  All  particle  properties  at  the  exit  plane 
ol  the  rocket  motor  including  particle  size  are  generated  by  SPP.  SPP 
provides  three  particle  size  groups  with  a  very  limited  amount  of  size  control 
in  the  hands  of  the  user.  Considerable  controversy  surrounds  the  actual  size 
distribution  of  the  particles  and  whether  Che  three  size  groups  generated  by 
SPP  do  or  do  not  give  an  adequate  model  of  the  particle  flow  field.  SCIPPY 
uses  the  output  of  SPP  to  generate  the  flow  field  within  the  diffuser  and  any 
shortcomings  of  SPP  in  terms  of  the  particle  flow  field  are  propagated 
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throughout  the  diffuser.  SCIPPY  in  turn  introduces  its  own  problems  with 
regard  to  the  particles.  If,  as  discussed  in  the  preceding  section,  the 

particles  are  assumed  not  to  adhere  L->  the  wall,  they  will  accumulate  as  a 
debris  layer  along  the  wall  or  they  will  be  reentrained  within  the  flow. 

Within  SCIPPY,  particles  that  strike  the  wall  are  allowed  to  pass  through  the 
wall  and  are  lost  from  the  flow.  In  addition,  SCIPPY  will  at  some  point  in 
the  flow  drop  a  Mach  disk  and  will  lose  all  particles  which  pass  through  the 
Mach  disk.  The  radiation  heat  load  is  a  minor  threat  to  the  diffuser  and 
until  the  particle  flow  field  is  better  defined  a  more  refined  model  does  not 
seem  justified. 

It  is  assumed  that  the  flow  field  is  optically  thin  and  that  the 
particles  behave  as  gray  bodies  emitting  radiation  as  per  the  Stefan-Boltzman 
equation  with  all  properties  evaluated  in  terms  of  centerline  conditions  at 
the  exit  plane  of  the  rocket  motor.  It  is  further  assumed  that  this  emissive 
power  is  concentrated  as  a  line  source  of  uniform  strength  along  the 
centerline  of  the  diffuser.  This  source  strength  is  readily  evaluated  in 
terms  of  exit  plane  information  from  SPP  and  takes  the  form 

q  =  3  m  coT^/u  OR  (5) 

p  p  p  p  p 
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wli?re  o  is  the  Stefan-Boltzman  constant  and  m  is  the  mass  flow  rate,  e  the 

P 

emissivity,  T^  the  temperature,  the  axial  velocity,  p  the  mass  density, 
and  Rp  the  particle  radius  associated  with  the  particle  group  in  question. 
Equation  5  must  be  summed  over  the  particle  groups  present.  Defense  of  this 
model  as  used  within  REDTAP  was  supported  by  the  assumption  of  a  thermal 
accommodation  coefficient  of  unity.  This  resulted  in  a  particle  impingement 
heat  Load  so  large  as  to  render  the  radiation  load  negligible.  With  the 
accommodation  coefficients  suggested  in  the  present  study,  the  radiation  load 
will  become  a  significant  but  not  major  portion  of  the  heat  load  on  the 
diffuser.  On  the  other  hand,  the  present  assumption  that  the  particles  do  not 
adhere  to  the  wall  lends  credence  to  the  assumption  of  a  uniform  strength  line 
source  of  radiation.  As  a  body  of  experimental  data  becomes  available  and  a 
revised  SCIPPV  provides  more  reliable  flow  field  data,  refinement  of  this 
model  should  be  considered. 
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?,4  GAS  SIDE  CGiJVSCTIOK 


Ine  gas  side  <onvecticn  heat  transfer  is  handled  within  DHT  by 
i  n«.  or  porn  t  in>;  TRL  as  a  subroutine  in  much  the  same  fashion  as  was  done  in 
KEDTAP.  SCIi’i’Y  provides  TBL  with  the  edge  conditions  for  the  boundary  layer 
analysis  and  TBL  provides  DHT  with  the  film  coefficients  and  adiabatic  wall 
temperatures  required  f.'r  the  heat  transfer  calculation,  The  boundary  layer 
erows  from  a  stagnation  region  at  the  poi.-.t  of  plume  impingement  and  therefore 
it  is  necessary  to  start  the  boundcry  Jayer  with  nonzero  initial  boundary 
layer  thicknesses.  REDTAP  started  TBL  with  the  initial  momentum  and  energy 

thicknesses  set  equal  to  a  single  arbitrarily  small  number.  This  approach  may 
provide  meaningful  inf 'rmaticn  well  removed  from  the  point  of  plume 
impingement  but  there  is  no  reason  to  expect  it  to  provide  useful  film 
coefficients  in  the  vicinity  of  plume  impingement.  Unfortunately,  this  is 

probably  the  most  critical  region.  In  addition,  no  attempt  was  made  to  handle 
the  stagnation  point  neat  transfer  which  occurs  at  the  poirt  impingement.  For 
some  cases,  Shis  may  provide  the  severest  heat  load  on  the  diffuser.  A  more 
reasonable  approach  is  Co  evaluate  the  impingement  point  heat  transfer 
coefficient  independent iy  and  to  start  the  boundary  layer  thicknesses  with 
values  that  will  result  in  this  level  of  heat  transfer.  This  approach 
admittedly  irop^pos  £  boundary  layer  model  on  a  region  that  is  not  a  boundary 
layer  flow  and  one  must  question  the  validity  of  the  results  immediately 

dowrstream  of  impingement.  Preliminary  experimental  results  do,  however,  show 
surprisingly  good  correlation  with  the  predictions  of  DHT.  This  is  the 
approach  used  within  DHT  and  is  in  essence  the  approach  suggested  by 

Pergsraent^. 

8 

Pergament  cites  the  work  of  Donaldson,  Snedeker,  and  Margolis  who  present 
data  for  a  0.51 1-inch  diameter  converging  nozzle  driven  by  compressed  air, 
discharging  to  the  atmosphere,  operations  at  pressure  ratios  of  O.dOO  to 
0.14b,  and  impinging  normally  upon  a  flat  surface.  Donaldson,  et  al,  show 
correlation  between  their  data  and  the  laminar  heat  transfer  correlation 

h  =  k  sT?TuU/d7)7Zv-  (6) 

where  H  is  the  film  coefficient,  k  is  the  thermal  conductivity  of  the  gas,  Pr 
is  the  Prandtl  Number  of  the  gas,  dU/da  is  the  velocity  gradient  of  the  free 
stream  flow  parallel  to  the  surface,  and  v  la  the  kinematic  viscosity  of  the 


gas — all  evaluated  just  downstreim  of  impingement .  The  film  coefficient  as 
defined  by  Equation  b  is  multiplied  by  an  augmentation  factor  which  ranges 
from  1.4  to  2.2  and  which  Donaldson,  et  al,  show  to  be  a  l unction  of  the 
distance  from  the  nozaie  to  the  impingement  surface  and  attribute  to 
turbulence  growth  within  the  jet  prior  to  impingement. 

While  there  would  be  no  reason  to  expect  the  work  of  Donaldson,  et  al,  to 
apply  directly  to  the  diffuser,  there  would  seem  to  be  hope  tot  using  Equation 
6  in  conjunction  with  an  augmentation  facto?  evolved  from  experimental 
diffuser  data.  The  velocity  grauitnt  dU/de  can  oe  evaluated  from  velocity 
data  supplied  by  SCIPPY.  While  this  is  an  area  of  the  tin-**  that  SCIPPY 
handles  ratner  poorly,  the  velocity  gradient  enters  as  a  square  root  and  the 
associated  er»-oi  is  reduced.  Preliminary  experimental  data  from  two  Super 
HATES  firings  show  good  agreement  with  predictions  mads  using  an  augmentation 
factor  of  2.5,  however,  an  understanding  of  this  area  will  come  only  as  we 
acquire  a  base  of  experimental  diffuser  data. 
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2.5.  DEBRIS  LAYER  SHIELDING 

The  assumption  that  the  particles  which  strike  the  wall  do  not  adhere  to 
the  wall  gives  rice  to  an  accumulation  of  these  particles  in  the  vicinity  of 
the  wall.  As  this  accumulation  is  swept  downstream  by  the  main  flow,  i.t  will 
form  an  increasingly  dense  sheath  of  particles  adjacent  to  the  wall  and  will 
partially  shield  the  wall  from  particle  impingement.  Wickoan,  Mockenhaupt , 
and  Ditore  develop  a  simple  model  for  this  phenomenon  and  present  supporting 
data  in  conjunction  with  an  erosion  study.  The  essentials  of  the  model  are 
contained  in  Figure  2.  The  model  assumes  a  single  particle  sire  and  a  cross 
sectionaL  area  for  collision  equal  to  o  .  Assuming  a  particio  number  density 
n  within  the  debris  layer,  the  cross-sectional  area  blocked  p**r  unit  area  by 
the  debris  is  found  to  be 
n"dx  /  SIN  B. 

Assuming  an  incident  particle  flow  with  a  particle  number  density  of  N,  the 
change  in  particle  number  density  caused  by  scattering  within  the  debris  layer 
element  dx  will  be 

dN  =  Nnodx  /  SIN  B. 
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Substituting  this  expression  into  Equation  ?  one  obtains 


Iihnw/n0; 


as./  Z  '  Rja  U  ,  SIN  B 
o  P  g 


It.  the  case  of  diffuser  flow,  md  is  obtained  by  summing  the  particle  mass 
flux  impinging  upstream  of  the  point  in  question.  Information  necessary  for 
evaluating  everything  except  o  is  available  from  SCIPPY. 

The  model  just  described  is  readily  expanded  to  include  flows  involving 
more  than  a  single  size  particle.  It  is  convenient  to  treat  each  particle 
size  group  individually.  Assuming  three  3ize  groups  with  group  j  assumed  to 
be  the  incident  particle  group  and  group  k  the  particle  debris  group,  one  can 
consider  as  the  cross  section  for  particle  group  j  colliding  with  particle 
group  k,  V.  the  particle  number  density  of  the  impinging  particles,  and  n^  the 
particle  numccr  density  of  the  debris.  With  this  nomenclature,  Equation  7 


may  be  written  as 


'"Wj-  L  -  / 


n^dx  / 


SIN  B. 


where  H  is  sufficiently  large  to  include  all  three  debris  layers.  The  righc 
hand  side  of  Equation  9  may  be  expanded  and,  noting  that  as  x  tends  to  d,  n 
tends  to  zero,  the  upper  limit  of  each  integral  may  be  replaced  with  the 
individual  debris  layer  thickness,  leading  to  the  following  form: 


,„(« y ....  r** 
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n^x/siN  Bj. 


As  with  Equation  S,  it  is  convenient  to  recognize  that 


o  j  \  ds  ^  \  1  2  KVg 

where  ta.^  is  the  mass  flow  rate  of  particle  group  k  within  the  debris  layer  and 
is  the  mass  of  a  group  k  particle.  Introducing  this  into  Equation  10  leads 
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Noting  that,  with  the  exception  of  SIN  B  j ,  the  above  expression  is  solely  a 

function  of  the  debris  layer,  it  is  convenient  to  define  DP,  the  debris 
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factor,  such  that 


ln( DFj )  = 


;  ( »jkvv  1 2 '"“e 


(HA) 


and 


<VV  =  (DFj ) 1/SIK  Bj 
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In  the 


No  mention  has  been  made  thus  far  as  to  evaluating  jk* 
development  of  their  particle  diameter  model,  Wickman,  et  al,  assume  chat  any 
contact  at  ail  between  impinging  particle  and  debris  particle  will  result  in 
the  scattering  of  the  impinging  particle.  This  model  leads  to 
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which  would  appear  to  be  excessive.  The  model  built  into  DHT  assumes  that  a 
smaller  particle  will  be  scattered  oy  as  little  as  grazing  contact,  that  an 
equal  size  particle  will  require  an  angle  of  impact  cf  at  least  45  degrees, 
and  that  a  larger  particle  must  impact  a  smaller  particle  with  an  angle  of  at 
least  45  degrees  and  impact  an  aggregate  mass  of  such  particles  equal  to  its 
own  mass  before  scattering  will  occur.  This  leads  to 
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Introducing  this  monel  for  o  ..  and  noting  that 
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one  is  able  to  evaluate  the  individual  terms  of  the  right  hand  side  of 
Equation  if A  as  f 

*  "  R.  <  B. 
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The  above  debris  shielding  model  has  been  built  into  DliT  and  the  user  is 
provided  with  the  option  to  use  or  not  use  it  in  the  calculations.  If  the 
option  is  implemented,  the  particle  mass  flow  rate  that  appears  in  Equation  3 
will  be  multiplied  by  the  factor 
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For  toe  examples  looked  at  to  date,  debris  shielding  has  not  appeared  to  be  a 
significant  factor.  The  debris  factor  has  ranged  from  1.0  to  0.9  but  nas 
remained  very  close  to  1.0  ir.  the  regions  where  impingement  heating  was  a 
major  concern.  This  is  understandable  since  only  after  particles  impinge  upon 
the  wall  for  some  distance  does  the  debris  layer  build  up  to  an  effective 
shield.  The  reduction  of  particle  mass  flux  teaching  the  wail  may  be  as  great 
as  30  percent  in  some  regions  but  these  regions  are  well  downstream  of  the 
severe  heat  load  areas.  The  regions  where  appreciable  debris  shielding  occurs 
are  where  the  impingement  angle  is  quite  shallow  and 
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can  become  quite  small.  For  sue4-  flows,  however,  if  one  is  assuming 
accommodation  coefficients  of  0,  0,  and  0.8  SIN  B,  the  impingement  heat  load 
is  quite  small  with  or  without  debris  shielding. 


2.6  WATER  SIDE  CONVECTION 

The  water  side  film  coefficient  is  evaluated  using  correlations  presented 
by  Marks.  *0  The  preliminary  calculation  is  handled  by 

h 1  -  160  (1  +  0.012  T,)  V°-W-2  (13) 

t  h 


i 


where  Tf  is  the  film  temperature  of  the  water  measured  in  degrees  Fahrenheit, 
V  vs  the  velocity  of  the  water  measured  m  feet  per  second,  is  the 
hydraulic  diameter  (A  x  area  /  perimeter)  of  the  channel  measured  in  inches, 
and  h*  is  the  film  coefficient  measured  in  B/hr-ft  -F.  This  value  for  h  is 
modified  to  compensate  for  the  radius  of  curvature  of  the  channel  such 
that 


h  »  (l  +  3.5  (D. /D  >)  h’ 
h  c 


(14) 
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3.  NUMERICAL  ANALYSIS 


Figure  3  shows  the  finite  difference  grid  that  is  used  in  solving  for  the 
temperature  distribution  within  the  diffuser  wall  and  the  water  jacket. 


where  the  superscript  +  indicates  a  temperature  occurring  at  time  (t  +  At), 


Equation  15  may  be  solved  fo»'  T  and  written  as 

rat  ti 


_  (T  .  +  T  .  +  Z (T  .  +  V  . .  ) 

i,n  =  m,n-l  m,n+l  m-l,n  m-l,n 


♦  (Ml  -  2  -  L  Z')  T  )  /  Ml 

m  t  n 


where 


2  =  i  *  Y /:  X )  CCS  0 


Ml  n  o  C(<-  Y)  *  /k  At . 


Ic  this  !orm  the  temperature  distribution  at  tisr-t;  (t  +  t)  may  be  solved  for 
point  by  point  in  terms  of  a  known  temperature  distribution  at  time  t.  This 
explicit  formulation  has  the  stability  requirement  thvt 


Ml  -  2  -  2  r  >  0 


For  a  given  AX  and  AY  this  places  an  upper  bound  on  At  but  has  presented  no 
prooletns  to  date. 

Equation  16  is  applicable  to  all  internal  nodes.  However,  the  first 
radial  node  (m,l)  and  the  lart  radial  node  (m,L)  involve  boundary  conditions 
and  must  be  handled  separately.  In  the  case  of  these  two  nodes  it  is 
convenient  to  forsake  the  mathematical  elegance  of  finite  difference  forms  and 
to  perform  an  energy  balance  on  the  element.  Kote  that  in  terms  of  thermal 
capacity  each  of  these  nodes  involves  only  one  half  an  element.  Written  in 
explicit  for  node  m,l  this  takes  on  the  following  form: 


h  Xl TAW  -  T  .)  ♦  UAY/2AX)  (T  ,  .  -  T  ,) 
m  m  m,l  m-1,1  m,l 


+  ( k Ay / 2 AX )  (T  -  T  .}  +  (kAX/AY)  (T  ,  -  T  ,) 

m+l,l  m,l  m,2  m,l 


+  QPRAX  (TP  -  T  .)  +  QPI  AX  +  QPR  AX 
m  m  m,l  m  m 


*  (pCAXAY/2 At) (T  ,  -  X  ,) 
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where  h  is  the  gas  side  film  coefficient,  TAW  is  the  adiabatic  wall 
temperature  of  the  gas,  QPT  is  (.he  thermal  energy  flux  per  unit  area  caused  by 
particle  impingement,  QP1  is  the  inertial  energy  flux  par  unit  ares  caused  by- 
particle  impingement,  and  QPR  is  the  radiant  energy  flux  per  unit  area  from 
the  particles.  Equation  17  may  be  solved  for  T+  ,  and  written  as 
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tii,  l  ra-l,l  m+ltl  m 
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tn  ra  id  m 


+  (Ml-  2  -  27.  -  2N1  -  2QPT  A  Y/k)  T  ,  )/Ml 
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(18) 


where 


hi  -  hAY/k. 


Here  as  with  Equation  lb,  one  has  the  simplicity  of  an  explicit  formulation 
but  the  stability  restriction  that 
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Ml  -  2  -  2Z2  -  2N1  -  2QPT  A  Y/k  >  0 
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which  places  a  slightly  smaller  upper  limit  on  A  t  than  was  associated  with 
Equation  16. 

Node  m,L  may  be  handled  in  tne  same  fdohion  as  node  m,l  and  results  in 


T+  ,  »  (2  T  ,  .  +  Z  (T  ,  ,  +  T  .  ,)  +  2N2TC 
m,L  m,L-l  ra-l,L  ra+l,L  m 
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where  TC  is  the  local  coolant  temperature  and  N2  is  identical  to  N1  only  based 
upon  the  water  siae  film  coefficient.  In  this  case  stability  requires  that 


Ml  -  2  -  2Z  -  2K2  >.  0 


r‘tich  places  an  additional  upper  bound  on  At 


SH 


An  attempt  to  handle  Equation  2  in  the  same  fashion  ss  Equation  1,  that 
is  to  say,  using  a  central  difference  approximation  for  the  spatial  derivative 
and  a  forward  difference  approximation  for  the  temporal  derivative  will  lead 
to  numerical  instability.  On  the  other  hand,  using  a  backward  difference 
approximation  fo.  the  spatial  derivative  will  lead  to  a  stable  formulation. 
Using  this  latter  approach,  Equation  2  may  be  approximated  by 

Um  <TCm  ~  TW/  ''X  +  (2"  Wp  CV(TCm  "  Tm,L> 


-  (TC«  -  TCm>/At 


Solving  for  TC 


where 


TC  =  ilC  ,  *•  N3  T  ,  +  (M2  -  I  -  N3)  TC  ]  /  M2 
ns  '  •  in,  L  in 


M2  »  A  X/U  A  t  =  oA  A  X  /  m  A  l  COS 
m  m 


N3  =  2sR  h  AX/roC  COS  0 
m  m 


and  n>  is  the  mass  flow  rate  of  the  coolant.  Stability  will  require  that 
M2  -  l  -  N3  >  0 

and  place  yet  another  upper  bound  on  Ac. 


4.  IMPLEMENTATION  OF  DHT 


The  following  discussion  is  intended  to  be  an  overview  of  the  code. 
Detailed  discussion  of  the  input  and  output  will  appear  separately. 

Prior  to  executing  DHT,  SSP  must  be  run  for  the  rocket  motor  and  SCIPPY 
must  be  run  for  the  diffuser.  In  running  SPP  it  is  necessary  to  save  and 
catalog  the  exit  plane  information  contained  on  TAPE29.  This  information  and 
a  small  amount  of  hand  input  data  is  required  by  SCIPPY.  In  running  SCIPPY  it 
is  necessary  to  catalog  ana  save  the  data  contained  in  TAPE99.  Data  contained 
in  this  file  includes  the  diffuser  geometry  along  with  the  edge  conditions 
required  for  the  evaluation  of  the  heat  loads  on  the  diffuser. 

Units  for  input  data  are  selected  on  the  basis  of  user  convenience  and  are 
converted  internally  to  pounds  force,  feat,  seconds,  degrees  Rankine,  and 
pounds  mass.  Units  for  output  are  selected  on  the  basis  of  user  convenience. 

The  edge  condition  data  provided  by  SCIPPY  on  TAPE99  are  randomly  spaced 
along  the  axis  of  the  diffuser.  Since  the  numerical  analysis  within  DHT 
assumes  a  uniformly  spaced  grid  system,  the  first  major  operation  within  DHT 
is  to  read  the  SCIPPY  tape  and  to  interpolate  within  the  data  to  create  a  set 
of  uniformly  spaced  edge  conditions.  In  conjunction  with  this  manipulation, 
all  edge  condition  type  calculations  within  the  DHT  model  are  also  performed. 
These  include  the  calculation  of  particle  impingement  mass  flux,  debris  layer 
data,  and  particle  related  heat  fluxes.  In  as  much  as  neither  the  rocket 
motor  chamber  pressure  nor  the  test  cell  pressure  remain  constant  throughout 
the  test,  provision  is  made  to  update  the  edge  conditions  by  way  of  additional 
SCIPPY  tapes.  This  provision  is  handled  through  the  index  SCIPPY  and  numbered 
data  sets  identified  as  SCIP(l)  through  SCIP(SCIPPY) . 

Once  the  edge  conditions  have  been  established,  DHT  prepares  to  call  TBL 
in  order  to  obtain  the  gas  side  film  coefficients  and  adiabatic  wall 
temperatures  of  the  gas.  The  input  data  has  included  a  value  for  the  film 
coefficient  at  the  point  of  plume  impigement,  obtained  from  Equation  6  and 
adjusted  by  an  augmentation  factor.  It  is  necessary  for  DHT  to  establish 
values  for  the  boundary  layer  thicknesses  at  the  point  of  plume  impingement 
such  that  TBL  will  predict  this  value  for  the  impingement  point  film 
coefficient.  DHT  will  accept  a  user  selected  initial  trial  value  for  the 


energy  thickness,  hold  the  momentum  thickness  equal  to  a  user  selected 
constant  time."  the  energy  thickness,  and  iterate  to  seek  a  set  of  momentum  and 

energy  thicknesses  that  will  result  in  the  desired  impingement  point  film 
coefficient.  At  best,  this  is  no  more  than  a  means  lo  an  end.  It  will, 

however,  result  in  an  impingement  point  heat  flux  in  keeping  with  a  known 
correlation.  There  is  no  reason  to  expect  this  to  predict  valid  film 
coefficient  data  immediately  downstream  of  impingement  but  tyopefully  local 
edge  conditions  will  dominate  and  produce  valid  results  further  downstream. 
The  validity  of  this  approach  will  have  to  be  reviewed  as  experimental  data 
become  available.  Since  the  film  coefficients  predicted  by  TBL  are  mild 
functions  of  the  gas  sid.*  surface  temperatures,  provision  is  made  through  the 
parameters  ICALL  and  DCALL  to  update  the  film  coefficients  as  the  temperature 
of  the  diffuser  wall  rises.  The  code  has  an  initial  update  built  in  to  it 
that  occurs  ten  time  increments  into  the  calculations.  Subsequent  to  this,  an 
update  occurs  every  DCALL  time  increments  or  with  each  new  set  of  edge 
conditions  obtained  by  way  of  a  SCIPPY  tape.  TBL  is  a  time  consuming  code  and 
indiscriminate  updating  should  be  avoided.  To  date  it  has  been  adequate  to 
update  TBL  just  prior  to  the  final  output  and  verify  that  no  significant 
changes  have  occured  in  the  film  coefficients. 

Many  of  the  coefficients  within  the  governing  equations  are  independent  of 
temperature  and  are  evaluated  as  preliminary  calculations.  Much  of  this 
information  is  output  as  a  matter  of  convenience  to  the  user. 

The  output  that  occurs  at  time  zerr  or  following  the  reading  of  a  new 
SCIPPY  tape  is  quite  extensive  and  contains  a  great  deal  of  boundary  layer, 
debris  layer,  and  individual  heat  flux  information  that  will  remain  constant 
throughout  the  calculations.  Future  output  will  occur  every  DOUT  time 
increments  and  is  appreciably  more  abbreviated. 

If  the  parameter  DAFLAG  has  been  read  in  as  other  than  zero,  DHT  has  the 
capability  to  store  gas  side  wall  temperature,  water  side  wall  temperature, 
coolant  temperature  as  functions  of  either  time  or  position  for  future  use  in 
creating  graphic  information.  Unless  this  option  is  implemented,  temperatures 
are  stored  for  time  t  and  (t  +  At)  only.  A  more  extensive  discussion  of  the 
option  is  contained  under  the  separate  heading  of  GRAPHIC  OUTPUT. 

The  main  heat  transfer  calculation  is  an  implementation  of  the  finite 
difference  equations  presented  earlier  in  Section  3.  It  is  in  order  to 


comment  on  the  handling  of  the  end  conditions.  The  near  end  and  far  end  of 
the  diffuser  wall  are  treated  as  adiabatic  planes.  These  conditions  are 
implemented  by  extending  the  grid  system  one  grid  line  beyond  each  end  and 
step  by  step  assigning  the  outboard  nodes  mirror  image  values  from  internal 
nodes.  This  allows  nodes  on  the  two  end  planes  to  be  ha  died  as  though  they 
were  internal  nodes  and  introduces  no  additional  equations.  The  coolant 
temperature  at  the  inlet  plane  of  the  diffuser  is  held  at  the  supply 
temperature  and  the  upstream  differencing  used  within  Equation  20  requires  no 
knowledge  of  coolant  temperature  beyond  the  exit  plane  of  the  diffuser. 

Following  each  time  step  there  is  the  opportunity  to  output  the 
temperature  distribution  through  the  parameter  DOUT,  to  store  a  portion  of  the 
temperature  data  under  the  DAFLAG  option,  to  update  the  gas  side  film 
coefficients  u^dcr  the  DCALL  option,  and  to  update  the  edge  conditions  with  a 
new  SGIPPY  tape.  Independent  of  any  of  the  above  mentioned  options,  the  water 
side  film  coefficients  are  a  function  of  the  average  film  temperatures  and  are 
updated  following  each  time  step. 


5.  GRAPHIC  OUTPUT 

Provision  is  made  for  saving  data  seta  for  use  with  the  Tektronirs  EZGRAF 
graphic i  system.  This  option  is  implemented  by  netting  the  flag  DAFLAG  to  any 
non.-.ero  integer.  Data  sets  will  be  written  to  TAPE99.  The  first  data  set 
written  will  be  a  set  of  axial  locations  measured  in  inches  and  locating  grid 
lines  on  an  index  interval  of  DXBA.  Following  this  data  Set  will  be  a  series 
of  three  data  sets,  oach  consisting  of  a  single  time  followed  by  a  set  of  gas 
side  wall  temperatures,  water  side  wall  temperatures,  or  coolant  temperatures 
occurring  at  axial  locations  consistent  with  the  initial  set  of  axial  location 
data.  A  series  of  three  such  data  sets  will  be  written  on  the  same  time 
interval  as  the  printed  output  from  the  code.  The  next  data  set  written  will 
be  a  time  base  for  displaying  temperature  as  a  function  of  time.  A  sequence 
of  times  will  be  written  on  a  time  interval  of  DDAOUT  time  steps.  This  time 
base  is  followed  by  a  series  of  data  sets  consisting  of  a  single  axial 
location  followed  gas  side  wall  temperatures,  water  side  wall  temperatures,  or 
coolant  tempera  :ures  on  the  above  mentioned  time  interval.  As  many  as  five 
axial  locations  for  such  temperature  sets  may  be  selected  by  specifying  the 
axial  grid  index  MDATA-  All  data  sets,  with  the  exception  of  axial  position 
and  time,  are  written  m  rounded  integer  form  because  of  the  real  number 
storage  limitations  within  EZGRAPF.  Axial  position  and  time  are  retained  in 
real  number  form. 
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6.1  Nomenclature 


Variable 

ACN 

ACP 

ACT 

ATHETA 

DAFLAG 

DCALL 

DDAOUT 

D1SCH 

DHXBA 


% 

X  1 

DOUT 

DX 

> 

I 

t 

DXMAX 

1 

! 

-  -  4 

DTAU 

r’  -  * 

* 

•  j  - 

•i  - 
1  : 

ENDTBL(K) 

—  1 1 

*■  •iL-j  1 

GAMO 

*  1  1 
’■  -  te  a  % 

HIMPK(K) 

Description 

Accommodation  coefficient,  kinetic  energy  normal 
to  the  diffuser  wall  (none) 

Accommodation  coefficient,  kinetic  energy  parallel 
to  the  diffuser  wall  (none) 

Accommodation  coefficient,  thermal  energy  (none) 

Proportionality  constant  (none) 

THETA I  =  ATHETA  *  PHI I 

Flag.  Data  sets  for  graphic  output  are  stored 
if  DAFLAG  is  nonzero,  (none) 

Frequency  of  TBL  update,  every  DCALL  time  steps  (none) 

Frequency  with  which  temperature  is  saved  for  graphic 
output  as  a  function  of  time,  every  DDAOUT  time  steps 
(none) 

Volumetric  flow  rate  of  coolant  (gpm) 

Frequency  with  which  temperature  is  saved  for  graphic 
output  as  a  function  of  time,  every  DHXBA  grid  lines 
(none) 

Frequency  of  printed  output,  every  DOUT  time  steps 
(none) 

Axial  seep  size  (inches) 

Maximum  allowable  axial  step  size  within  TBL 
( inches) 

Time  step  size  (sec) 

Last  time  step  for  which  the  Kth  SCIPPY  tape  should 
be  used  (none) 

Stagnation  ratio  of  specific  heats  associated  with 
the  Kth  SCIPPY  tape  (none) 

Gas  phase  heat  transfer  coefficient  at  the  point 
of  plume  attachment  associated  the  Kth  SCIPPY  tape 
(B/sec-ft2-R) 


Type 

Reai 

Real 

Real 

Real 

Integer 

Integer 

Integer 

Real 

Integer 

Integer 

Real 


I  i  ' 

?  » 

I  < 


T££C 

Real 


Var  iab  le 
HT 

KRPIK(K) 

XW 

MDATA(N) 

MU 

NCli 

NDAIA 

NDTAil 

NDY 

NSCIP 
PH  11 

PRK(X) 

RBARK(K) 

RHOC 

RHOW 

SC1?(K) 

SPHTC 

SPHTW 


Description 

Radial  height  of  the  water  jacket  passage  (inches) 


Radiation  source  strength  associated  with  the 

Kth  SCIPPY  tape  (B/sec-ft)  Real 

Thermal  conductivity  of  the  diffuser  wall 

(B/sec-ft-R)  Real 

Axial  grid  location  at  which  temperature  is  to  be 

saved  for  graphic  output  (none)  Integer 

Dynamic  viscosity  of  Che  coolant  (lbm/ft-sec)  Real 

Number  of  coolant  channels  (none)  Integer 

Number  of  axial  positions  MDATA(N)  at  which 


temperatures  are  to  be  saved  for  graphic  output  (none)  Integer 


Number  of  time  steps  of  numerical  analysis  to  be 

performed  (none)  Integer 

Number  of  diffuser  wall  elements  taken  radially 

(none)  Integer 

Number  of  SCIPPY  tapes  to  be  read  (none)  Integer 

Initial  estimate  of  the  energy  thickness  at  the 

point  of  plume  attachment  (feet)  Real 

Stagnation  Prandtl  Number  associated  with  the  Kth 

SCIPPY  tape  (none)  Real 

Gas  constant  associated  with  the  Kth  SCIPPY  tape 
( ft-lbf/ibm-R)  Real 

3 

Mass  density  of  the  coolant  (lbir/ft  )  Real 

Mass  density  of  the  diffuser  wail  (lbm/ff*)  Real 


Identifier.  Various  SCIPPY  tapes  may  be  attached 
as  TAPEnn.  SCIP(K)  is  the  two  digit  identifier  nn. 

Provision  is  made  within  the  DHT  program  card  for 

TAPEil  and  TAPE12.  This  may  be  expanded  if  the  user 

desires  and  has  the  storage  space  available,  (none)  Integer 

Specific  heat  of  the  coolant  (B/Jbm-R)  Real 

Specific  heat  of  the  diffuser  wall  (B/ib»-R) 


Real 


Variab  le 

Description 

Type 

T1 

Initial  temperature  (°R) 

Real 

TKW 

Thickness  of  the  diffuser  inner  wall  (inches) 

Real 

TOK(K) 

Stagnation  temperature  associated  with  the  Kth  SCIPPY 
tape  (°R) 

Real 

TYPACH 

Type  AON 

TYPACN  t  0,  ACN  =  ACN 

TYPACN  =  0,  ACN  =  0.8  SIN  B  (none) 

Integer 

TYPDBR 

Type  dehris  layer  analysis 

TYPDBR  =  0,  effects  of  model  excluded 

TYPDBR  =  0,  effects  of  model  included 

Parameters  will  be  calculated  and  outputted 
regardless  of  the  value  of  TYPDBR  (none) 

Integer 

XMOTOR 

Distance  from  the  motor  exit  plane  to  the  diffuser 
inlet  plane.  If  the  exii"  '.one  of  the  motor  extends 
into  the  diffuser,  XMOTOR  will  be  negative,  (inches) 

Real 

XSTOP 

Extent  of  diffuser  to  be  analyzed  as  measured  from 
the  inlet  plane  of  the  diffuser  (inches) 

Real 

WIDTH 

Width  of  each  coolant  channel  (inches) 

Real 

ZMOUKfK) 

Stagnation  viscosity  associated  with  the  Kth 

SCIPPY  tape  ( Ibm/f t-sec) 

ZHVISK(K; 

Exponent  associated  with  the  Kth  SCIPPY  tape  viscosity 
vs  temperature  model 

ZMVISK  =  ZMOUK  (T/T0K)  (none) 

Integer 

6.2  INPUT  PROCEDURES 

DHT 

receives  its  input  in  the  form  of  a  series  of  data 

card 

avai  la'ole 

as  TAPE5.  Shown  below  is  a  card  by  card  description  of 

the 

format 

and  data 

contained  within  each  card  image. 

CARD  1  U2A6) 
t ITLE  FOR  THE  ANALYSIS 

CARD  2  (8P10.4) 

DX,  DTAU,  DXMAX 


CARD  3  18110) 

NDTAU,  NDY,  NCH,  DOUT,  DCALL,  NSC IP 


I 


CARD  4  (8110) 

DAFLAG,  DDAOUT,  DMXEA,  NDATA ,  MDATA( 1 ) ,  - ,  MDATA(NDATA) 

If  no  data  sets  are  to  be  stored  for  graphic  output,  this  card  image  may 
be  left  blank  but  must  be  included. 


CARD  5  (8F10-4) 

TKW,  WIDTH.  HT 

CARD  b  ^SF!0.4i 

RHOW,  RHOC,  SPHTW,  KW,  HU 

CARD  7  (8F10.4) 

DISCH,  XHOTOR,  XSTOP,  TI 

CARD  8  (8F10.4) 

ACN.  ACP,  ACT 

CARD  9  (8110) 

TYPACN,  TYPDBR 

CARD  10  (8110) 

SCIP(i),  - ,  GCIP(NSCIP) 

CARD  11  (8110) 

ENDTBL(l),  - ,  ENDT3L(NSCl?) 

CARD  12  (8F10.4) 
l'OK(l),  - ,  TOK(NSCIP) 

CARD  13  (8F10.A/ 

ftBARX(l),  - ,  RBAKK ( NSC IP) 

CARD  14  (8F10.&) 

?RK( 1 ) ,  - ,  PRK(NSCIP) 

CARD  15  (8F10.4) 

ZMUOK(l),  - ,  2HUOR(NSCIP) 


CARD  16  (Sr  10.4) 

ZMVISK(i),  - ,  ZMVISK(SCIP) 

CARD  17  (8F10.4) 

C-AHOK(l), - ,  GAMOK(NSCIP) 


CARD  i8  (8F10.4) 

KRPIK(l),  - ,  KRPIK(NSCIP) 

CARD  19  (8F10.4) 

KIMPK(l),  - ,  HIMPK(NSCIP) 

CARD  20  (8F10.4) 

PHI I,  ATHETA 
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6.3  INPUT  GUIDELINES 

ATHETA,  PHI I. 

PHII  is  the  first  trial  value  of  the  energy  thickness  used  by  the 
iterative  procedure  that  seeks  initial  energy  and  momentum  thicknesses 
consistent  with  the  desired  impingement  point  heat  transfer  coefficient.  A 
value  of  PH11  that  is  either  too  large  or  too  small  will  cause  TBL  to  fail. 
PHII  equal  to  0.005  feet  has  been  found  to  be  successful  with  the  motors 
considered  to  date.  ATHETA  is  the  ratio  of  the  momentum  thickness  to  the 
energy  'hickness  and  is  used  only  in  the  iterative  procedure.  A  value  of  1.1 
has  been  found  to  be  successful  with  the  motors  considered  to  date.  A  value 
of  1.0  has  often  caused  TBL  to  fail. 

NSCIP,  SCIP  (K). 

For  simple  analyses  the  flow  field  within  the  diffuser  will  be  assumed  to 
remain  constant  with  respect  to  time  and  NSCIP  will  be  1.  If,  on  the  other 
nand,  either  the  chamber  pressure  of  the  motor  or  the  test  cell  pressure 
change  significantly  during  burn  time,  it  may  be  desirable  to  account  for  the 
resulting  changes  in  the  diffuser  flow  field.  This  would  require  the  use  of 
two  or  more  SPP  and  SCIPPY  runs  and  would  generate  two  or  more  SCIPPY  tapes. 
If,  for  example,  the  flow  conditions  were  quite  similar  during  the  early  and 
late  portions  of  the  burn  time  but  differed  significantly  during  the  middle  of 
the  burn  time,  one  could  consider  using  two  SCIPPY  runs.  NSCIP  would  be  read 
in  as  3.  The  two  SCIPPY  tapes  could  be  attached  as  TAPE11  and  TAPE12. 
SCIP(l)  would  be  read  in  as  11,  SCIP(2)  would  be  read  in  as  12,  and  SCIP(3) 
would  be  read  in  as  11.  ENTBL(i),  ENTBL12),  and  ENDTBL(3)  would  indicate  the 
last  time  step  for  which  each  set  of  data  would  be  used.  Cards  12  through  19 
will  contain  NSCIP  entries.  If,  as  in  the  example  just  cited,  a  given  tape  is 
to  be  used  more  than  one  time,  there  will  be  duplication  among  the  entries  but 
there  will  be  NSCIP  entries  per  card.  Provision  has  been  made  in  the  program 
card  for  TAPE11  and  TAPE12.  The  user  may  expand  on  this  as  machine  time  and 
space  permit.  The  program  is  dimensioned  to  allow  NSCIP  to  be  as  large  as  10. 

HIMPKUO,  XRPIKUO. 

See  sample  problems  for  details  on  calculating  these  parameters. 
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Variable 
PHI  I 


Description 


QPIJ(J) 


SCIPPY 


SINEJ(J) 


Energy  boundary  layer  thickness  at  the  point  of  plume 
impingement  (ft) 

2 

Water  side  convective  heat  flux  (B/sec-ft  ) 

2 

Gas  side  convective  heat  flux  (B/sec-ft  ) 

Inertial  heat  flux  associated  with  particle  impingement 
(B/sec-f t2) 

Inertial  heat  flux  associated  with  particle  impingement, 
group  J  only  (B/sec-ft^) 

Heat  flux  associated  with  particle  radiation  (B/sec-ft^) 

Thermal  heat  flux  associated  with  particle  impingement 
(B/sec-f c  2 ) 

Reynolds  Number  for  the  coolant  flow  (none) 

Local  diffuser  radius  (feet) 

Local  inner  radius  of  the  water  jacket  (feet) 

Local  outer  radius  of  the  water  jacket  (feet) 

Tha  number  of  the  SCIPPY  tape  being  used  (none) 

Sine  of  the  impingement  angle  with  which  particle  group  J 
strikes  the  diffuser  wall  (none) 


TAW  Adiabatic  wall  temperature  of  the  edge  condition  gas  flow 

(°P.) 

V  Coolant  velocity  (fps) 

X  Axial  location  measured  from  the  inlet  plane  of  the 

diwuser  (inches) 

Y  Diffuser  inside  radius  (feet) 

7.2  GENERAL  DESCRIPTION 

The  output  from  DHT  is  well  labelled  and  with  the  above  nomenclature 
should  be  self-explanatory;  however,  a  general  description  of  the  output 
should  be  useful  to  the  first  time  user. 

The  initial  set  of  information  generated  by  DHT  consists  of  a  user 
defined  title  for  the  analysis  followed  bv  a  listing  of  the  input  data 
supplied  by  the  user. 
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At  this  point  DIIT  will  read  and  organize  the  data  contained  on  the  SCIPPY 
tape  and  will  then  call  TBL.  The  next  several  sets  of  output  information  froc\ 
DHT  will  actually  be  generated  from  within  TBL.  The  first  of  these  sets  of 
information  will  be  a  set  of  boundary  layer  parameters.  This  will  be  followed 
by  a  set  of  geometry  data  (radius  vs  axial  location)  for  the  diffuser, 
freestream  Mach  number,  and  gas  side  wall  temperature  for  the  diffuser. 
Following  this  will  be  a  set  of  freestream  pressure,  temperature,  velocity, 
and  density  data  for  use  in  the  boundary  layer  analysis.  At  phis  point  the 
listing  will  only  show  data  for  the  first  four  axial  stations  following  the 
point  of  plume  impingement.  This  is  caused  by  the  fact  that  DKT  is  at  this 
point  setting  up  to  begin  an  iterative  routine  for  selecting  starting  values 
for  the  momentum  and  energy  tioundary  layer  thickness  that  will  result  in  TBL 
matching  a  user  supplied  heat  transfer  coefficient  at  the  point  of  plume 
impingement.  DHT  will  at  this  point  output  a  series  of  sets  of  PHII,  HG(1), 
and  MIMPK  tnat  will  allow  the  user  to  observe  the  process  by  which  DKT  will 
select  the  starting  boundary  layer  thicknesses.  Once  this  process  is 
complete,  DHT  will  call  TBL  and  generate  boundary  layer  information  for  the 
entire  diffuser.  This  will  generate  as  output  a  aw  seL  of  boundary  layer 
parameters  that  will  include  the  selected  starting  boundary  layer  thicknesses 
followed  by  a  complete  listing  of  diffuser  wall  geometry  and  edge  conditions. 

Ai  this  point  DHT  will  perform  a  number  of  preliminary  calculations  and 
output  preliminary  data  that  will  include  Ml,  RE,  V,  HC1,  Cl,  SCIPPY,  KRPIK, 
Dl'(K),  R(M),  R 1 ,  R2,  A,  M2,  and  C2.  During  the  listing  of  this  information, 
DHT  will  perform  a  check  on  Che  various  stability  criteria  and  print  a  warning 
message  if  any  are  being  violated. 

DHT  is  now  ready  to  begin  the  mam  heat  transfer  analysis.  The  initial 
condition  listing  that  appears  at  this  point  is  quite  comprehensive  and 
includes  numerous  parameters  that  will  remain  constant  thrcughout  the  entire 
analysis  or  at  least  until  a  new  SCIPPY  tape  is  read  in.  Information  will  be 
listed  at  every  axial  station  and  will  include  X,  MDOT,  SINE,  DEBRIS,  DEBRIS 
FACTOR,  MDOTW ,  QPIJ,  DELTA,  HG,  QHG,  QPI,  QPT,  QPR,  QHC,  TAW,  T,  and  TC.  It 
should  be  noted  that  the  wall  temperatures  are  labelled  as  WALL  and  are  listed 
sequentially  from  the  gas  side  wall  temperature  to  the  water  side  wall 
temperature.  Following  this  initial  listing  the  frequency  of  output  is 
controlled  by  the  parameter  DOUT  and  is  appreciably  abbreviated  unless  a  new 
SCIPPY  tape  is  read  in  which  case  the  more  comprehensive  listing  is  triggered. 


8.  SAMPLE  PROBLEM 


8 . 1  GENERAL 

DHT  lias  been  used  to  analyze  two  Super  BATES  firings  conducted  at  the 
AFRPL  Test  Area  1-42  and  to  provide  predictions  for  a  proposed  Minutetnan  III 
Stage  3  firing  11  cb?  same  f  •>'-?' ity.  The  necessary  input  data  for  the  January 
13,  1383  Super  3A"ES  firm,  will  be  present  '<  her*  along  with  the  DHT 

predictions  lor  all  three  firings 

8.2  Prc  i  j‘\ nary  Calculations 

Mt  o£  the  input  mfortroiicn  required  by  DHT  is  available  from  the 
«tatorsent  of  the  problem  or  the  output  from  3PP  or  SCIPPY.  There  are, 
however,  two  item.*  that  are  left  as  preliminary  hand  calculations.  The 
impingement  point  heat  transfer  coefficient  must  be  evaluated  from  Equation  6 
with  the  free  stream  velocity  gradient  approximated  from  the  edge  condition 
information  coi  tained  m  the  SCiPPY  tups. 

h  =  k  -\f  Pr(Cl)'/  ds)  /  2\> 

^  =  4.58  x  10"5  B/ sec-f t-°R 

Pi  «  0.4547 

dO/da  *  1006  sec 

v  “  0.0452  ft2/sec 

h  “  U. 00326  B/sec-ft2-°R 

Equation  6  is  a  laminar  correlation  and  must  be  tnltiplied  by  an  augmentation 
factor  which  for  the  moment  remains  unknown.  it  is  hoped  that  as  experimental 
diffuser  data  become  availaole  this  factor  will  be  *  time  better  defined  but  for 
the  moment  a  figure  of  2.5  is  viewed  as  a  conservative  estimate  and  leads  to: 

H1MPK  =  0.0082  B/sec-ft2-°R 

The  radiation  source  strength  KRPI  must  be  evaluated  using  Equation  5.  The 
information  necessary  to  make  this  calculation  is  obtained  from  the  output  of 
SPP.  SPP  partitions  the  mass  flow  of  particles  into  three  parcicle  groups 
based  upon  diameter.  These  three  particle  groups  each  carry  the  same  mass 
flow  rate.  Given  the  uncertainty  surrounding  the  particle  flow  field,  the 
simplistic  radiation  model  being  used,  and  the  relative  aic.e  of  the  radiation 
heat  flux,  it  seems  appropriate  to  evaluate  this  source  strength  for  the 


middle  size  pai  tide  group  and  to  multiply  this  figure  by  three.  At  such  a 
time  as  the  particle  flow  is  better  understood  it  would  be  appropriate  to 
revise  not  only  this  procedure  but  the  radiation  model. 


i 


q  =  )m  eel  /U  P  R 
P  P  P  P  P 


q  =  88  B/sec-ft 


=  16.72  lbra/sec 


8556  fps 


=  4139  R 


248  lbm/ft' 


=>  0.1714  x  10~8  B/hr-ft2-V> 


=  9.829  x  10  ft 


This  figure  is  the  source  strength  resulting  from  a  single  particle  group  and 
must  be  multiplied  by  three  such  that 
KRPl  a  264  B/sec-ft. 

8.3.  INPUT  DATA 

The  following  data  is  necessary  m  order  to  run  the  program. 


=  1.0  in 


=  0.05  sec 


OXMAX 


=  1.0  in 


NDTAU 


=  100  (this  will  provide  5  seconds  of  data) 


DCALL 


=  20  (this  will  provide  da?«  every  second) 

*  90  (rhia  will  update  TBi.  just  prior  to  the  end  of  the  run) 
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NSC1P 

21 

1 

(this  provides  for  using  only  one  SCIPPY 
tape) 

DAFLAG 

“ 

1 

(set  equal  to  non-zero  this  will  provide  for 
the  'ing  of  data  sets  for  use  with 
graphic  output) 

DDAOUT 

20 

(data  will  be  saved  every  second) 

DMXBA 

« 

1 

(data  will  be  saved  every  inch) 

NDATA 

S 

3 

(data  will  be  saved  as  a  function  of  time  at 
five  axial  locations) 

MDATA(N) 

= 

22,  23,  24,  25,  26 

(data  will  be  stored  as  a  function  of  time 
at  20,  21,  22,  23,  and  24  incne3) 

TKW 

= 

0.3  in 

WIDTH 

= 

5.25  in 

HT 

as 

2.75  in 

RHOW 

= 

490  lbm/ft3 

RHOC 

= 

62.4  lbm/ft3 

SPHTW 

5* 

0.1  B/  lbm-°R 

SPHTC 

C 

1.0  B/lbtn-°R 

KW 

= 

0.00863  B/sec-f t-°R 

MU 

= 

0.000759  lbm/ft-sec 

DISCH 

= 

1100  gpm 

XMOTOP, 

E 

12  in 

(the  exit  plane  of  the  motor  is  positioned 
12  inches  in  front  ol  tin-  inlet  to  the 
diffuser) 

XSTOP 

52  in 

(the  calculation  must  be  terminated  at  52 
inches  because  SCIPPY  fails  at  this  point 
for  a  motor  with  an  exit  plane  diameter  as 
small  as  Super  BATES  in  a  diffuser  as  large 
as  this) 

TI 

= 

510°R 

ACN 

= 

1.0 
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ACP  =  0.0 

ACT  =0.0 

TYI’ACN  =  0 

TYPDBR  =  0 

SC1P(N)  =11  (if  more  than  one  SCIPPY  tape  is  to  be 

called  this  will  be  a  sequence  of  tape 
numbers) 

ENDTBL(N)  =  100 

TOK(N)  =  6666°R 

ROARK (N )  =  30.13  i t- lb/  1  l\n~°R 

PkK(N)  =  0.4547 

ZMUOK(N)  =  0.00006741  lbm/f t-sec  * 

ZMVlSK(h)  =  0.656 

i 

GAMOK(N)  “1.29  * 

4 

KRPIK(N)  =  264  B/sec-it 

HIMPK(N)  =  0 . 0032  8/sec-f t-°R 

PH 1 1  =  0.005  ft 

ATliETA  =1.1 

8.4  EXECUTION  OF  THE  PROGRAM 

At  the  tune  of  execution,  the  above  data  must  be  available  to  DHT  as 
TAPES  and  in  a  format  as  specified  in  Section  6.  Listed  is  a  set  of  such  card 
images  consistent  with  the  above  data. 

CARD  1 

SUPER  BATES  -  13JAN83 
CARD  2 

1.0  0.05  1.0 

CARD  3 

100  4  4  20  90  1 

CARD  4 

1  20  l  5  22  23  24  23  26 
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(sSOS 
-  KSPiS 


sp 

ISSitest 


CARD  5 

0,5  5.25  2.75 


% 


! 


i  * 

:  t  ' 


CARD  6 


490.0 

62.4 

0.' 

1.0  0.00863 

0.000759 

CARD  7 
1100.0 

12.0 

52.0 

510.0 

CARD  8 

1 .0 

0.0 

0.0 

CARD  9 
0  0 


CARD  10 
11 


CARD  11 
100 

CARD  12 

6666.0 

CARD  13 
80.13 

CARD  14 
0.4547 

CARD  15 
0.00006241 

CARD  16 
0.656 

CARD  17 

1.20 


CARD  18 
264.0 

CARD  19 
0.0082 

CARD  20 
0.005  l.l 
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8 . 5  RESULTS 

Output  generated  by  DHT  and  based  on  these  data  is  presented  in  Figures  4 
and  5.  The  77-inch  diffuser  located  in  Test  Area  1-42  is  instrumented  to 
record  water  side  wall  temperatures  on  roughly  2-inch  centers  for  the  first  6 
feet  of  the  diffuser.  The  diffuser  is  also  instrumented  to  record  coolant 
temperatures  on  ro-  ;hly  the  same  intervals,  but  the  burn  time  of  the  Super 
BATES  motor  is  only  about  5  seconds  and  does  not  result  in  an  appreciable  rise 
in  the  coolant  temperature.  It  i«  unfortunate  that  SCIPPY  fails  after  the 
first  52  inches  of  the  diffuser  but  up  until  that  point  it  would  appear  that 
the  output  from  DHT  tracks  the  experimental  data  quite  well. 

SUPER  BATES  I3JAN83  -  TIME  =4.0  SECONDS 


520  j. 


-  THEORY 

-  INITIAL 

A  EXP  DATA  UfS 
v  EXP  DATA  RMS 


* . A  V  V  A  V  AV  A  9  A 
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Fjujure  4.  Axial  Temperature  Variations. 
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Figure  4  shows  the  axial  temperature  variations  of  the  water  side  wall 
temperature  at  approximately  the  end  of  the  burn  time.  Unfortunately  SCIPPY 
breaks  down  prior  to  the  point  at  which  there  is  any  particle  impingement  upon 
the  diffuser  wall.  Therefore  the  correlation  between  DHT  and  the  experimental 
data  seen  in  Fig.  4  in  no  way  validates  the  manner  in  which  the  code  is 
handling  particle  interactions  with  the  wall.  It  does,  however,  lend  support 
to  the  manner  in  which  the  code  is  handling  the  point  of  plume  impingement. 


SUPER  BATES  13JAN83 


Figure  5  shows  the  water  aide  wall  temperature  as  a  function  of  time  in  the 
vicinity  of  plume  impingement.  Here  again  the  degree  with  which  the 
predictions  of  OUT  track  the  experimental  data  is  very  gratifying,  but  the 
burn  time  of  the  Super  BATES  motor  is  insufficient  to  show  whether  the 
predictions  and  the  experimental  dati  will  compare  favorably  at  burn  times 
approaching  60  seconds.  Until  data  are  available  for  a  motor  with  a  burn  time 
approaching  at  least  30  seconds,  it  will  be  impossible  to  do  much  in  terms  of 
validating  the  predictions  of  DHT.  Also,  until  a  motor  which  SCIPPY  can 
handle  is  tested  or  until  a  better  version  of  SCIPPY  is  available,  very  little 
can  be  said  as  to  how  well  DHT  handles  the  particle  wall  interactions. 

Several  comments  are  in  order  concerning  the  experimental  data  shown  in 
Figure  A.  The  diffuser  is  instrumented  to  record  water  side  wall  temperature 
by  means  of  thermocouples  spot  welded  to  the  inner  wall  of  the  water  jacket. 
The  junctions  are  formed  by  spot  welding  each  thermocouple  lead  to  the  wall 
individually  and  allowing  the  wall  to  become  a  portion  of  the  thermocouple 
circuit.  It  was  hoped  that  by  so  doing  it  would  be  possible  to  locate  the 
effective  junction  at  the  surface  of  the  diffuser  wall.  It  is  felt  that  the 
lita  seer,  m  Fig.  4  stand  in  testimony  to  the  success  of  this  endeavor.  These 
thermocouples  are  sited  along  tw^  straight  line  paths  along  the  axis  of  the 
nffjser.  One  of  these  paths  labeled  LHS  is  located  roughly  45  degrees  from 
tiu*  top  and  along  the  left  hand  side  of  the  diffuser.  The  other  path  is 
roughly  45  degrees  from  the  top  and  along  the  right  hand  side  of  the  diffuser. 
Tne  data  shown  in  Figure  4  reveals  a  definite  biasing  of  the  data  in  terras  of 
’.eft  vs  fignt  hand  side  but  snows  either  set  of  data  from  a  single  side  of  the 
diffuser  to  be  very  self-consistent.  This  consistency  among  data  from  a 
single  side  of  'tie  diffuser  is  tele  to  rule  out  experimental  scatter  and  to 
validate  the  quality  of  the  data.  The  bias  seen  between  the  LHS  and  the  RHS 
is  definitely  real.  Th''  exact  same  trend  can  be  seen  in  data  gathered  from  a 
December  3,  1982  test  of  a  Super  BATES  motor  and  presented  in  Figure  1.  This 
same  trend  was  also  seen  in  an  entirely  different  motor  tested  in  this 
facility  prior  to  the  December  Super  BATES  firing.  This  left  vs  right 
variation  in  the  data  could  be  explained  in  terms  of  a  lack  of  symmetry  within 
the  diffuser  or  the  water  jacket  or  a  misalignment  of  the  test  stand  with 
respect  to  the  diffuser.  For  the  moment,  however,  it  is  real  and  measurements 
taken  on  the  facililty  reveal  no  misalignment  of  the  test  stand  with  respect 
to  the  diffuser. 


To  the  extent  possible  at  this  time  the  data  contained  in  Figs.  4  through 
6  support  the  validity  o£  both  the  code  and  the  experimental  data,  however, 
there  is  a  strong  need  for  experimental  data  from  a  motor  compatible  with 
SCIPPY  and  having  a  burn  time  of  30  seconds  or  longer. 

Figures  7  through  9  show  predictions  from  DHT  for  a  proposed  Minuteman 
III  Stage  3  test.  Unfortunately  this  test  has  yet  to  be  conducted  and 
therefore  there  are  no  experimental  data  for  comparison.  The  predictions  are 
included  as  being  indicative  of  the  capability  of  the  DHT  code. 
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TIME  »  0  SECONDS 
TIME  =  S  SECONDS 
TIME  -  10  SECONDS 
TIME  =  20  SECONDS 
TIME  °  60  SECONDS 
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